clear; clc; close all

dimw = 2;
dimx = 3;
dimy = 4;
dimz = 5;
n = 1000;
coeffs = [];

d(:,1) = round(rand(n,1));
d(:,2) = 5*rand(n,1) - 2.5;

for i = dimw+1:dimw+dimx,
    coeffs{i} = 5*rand(dimw,1) - 2.5;
    if(rand < 0.5)
        d(:,i) = d(:,1:dimw)*coeffs{i} + 2;
    else
        temp = exp(0.5 + d(:,1:dimw)*coeffs{i});
        d(:,i) = round(temp./(temp + 1));
        if(sum(isnan(d(:,i))) == 0)
            fprintf('%d:%d\n', i,sum(d(:,i)));
        else
            fprintf('ERROR!\n');
            return
        end
    end
end

for i = dimw+dimx+1 : dimw+dimx+dimy,
    coeffs{i} = 5*rand(dimw+dimx,1) - 2.5;
    if(rand < 0.5)
        d(:,i) = d(:,1:dimw+dimx)*coeffs{i} + 2;
    else
        temp = exp(0.5 + d(:,1:dimw+dimx)*coeffs{i});
        d(:,i) = round(temp./(temp + 1));
        if(sum(isnan(d(:,i))) == 0)
            fprintf('%d:%d\n', i,sum(d(:,i)));
        else
            fprintf('ERROR!\n');
            return
        end
    end
end

for i = dimw+dimx+dimy+1 : dimw+dimx+dimy+dimz,
    coeffs{i} = 5*rand(dimw+dimx+dimy,1) - 2.5;
    if(rand < 0.5)
        d(:,i) = d(:,1:dimw+dimx+dimy)*coeffs{i} + 2;
    else
        temp = exp(0.5 + d(:,1:dimw+dimx+dimy)*coeffs{i});
        d(:,i) = round(temp./(temp + 1));
        if(sum(isnan(d(:,i))) == 0)
            fprintf('%d:%d\n', i,sum(d(:,i)));
        else
            fprintf('ERROR!\n');
            return
        end
    end
end

d(1,find(d(1,:) == 0)) = 1;
d_orig = d;
d(1:n/2,dimw+1:dimw+dimx) = inf;
d(n/2 + 1:end,dimw+dimx+dimy+1:end) = inf;    

for i = dimw+1:dimw+dimx,
    randinds = find(rand(n/2,1) > 0.97) + n/2;
    d(randinds,i) = inf;
end

for i = dimw+dimx+dimy+1:dimw+dimx+dimy+dimz,
    randinds = find(rand(n/2,1) > 0.97);
    d(randinds,i) = inf;
end

% d(ceil(rand(round(numel(d)/50),1)*numel(d))) = inf;

d(:,1:dimw) = d_orig(:,1:dimw);

nMissingPerColumn = zeros(1,size(d,2));
for i = 1:size(d,2),
    nMissingPerColumn(i) = sum(d(:,i) == inf);
end

Ycoeffs = [1 2 3 4 0 0.1 0.2 -2 -8 -3 -2 -4 1.5 10]';
temp = exp(1 + d_orig*Ycoeffs);
Y = round(temp./(temp + 1));
sum(Y)
% return
data = [Y ones(n,1) d];
data_orig = [Y ones(n,1) d_orig];

dlmwrite('C:\Users\Mittal\Documents\Work\Data\Adolescent\Test\testAim2.txt',data,'delimiter',',');
dlmwrite('C:\Users\Mittal\Documents\Work\Data\Adolescent\Test\testAim2_orig.txt',data_orig,'delimiter',',');
csv2bbr('C:\Users\Mittal\Documents\Work\Data\Adolescent\Test\testAim2.txt');
csv2bbr('C:\Users\Mittal\Documents\Work\Data\Adolescent\Test\testAim2_orig.txt');
replaceinfile('Inf','NULL','C:\Users\Mittal\Documents\Work\Data\Adolescent\Test\testAim2.txt');
replaceinfile('Inf','NULL','C:\Users\Mittal\Documents\Work\Data\Adolescent\Test\testAim2_bbr.txt');
